function [dx dz x z ]=cont_grid(xmin,xmax,zmin,zmax,myx,myz,sigmax,sigmaz,I,H)

F=inline('0.5*(1+erf((x-my)/(sigma*sqrt(2))))', 'x','my','sigma');
Finv=inline('my+sigma*sqrt(2)*erfinv(2*x-1)', 'x','my','sigma');

x=zeros(1,I);
Fmin=F(xmin,myx,sigmax); Fmax=F(xmax,myx,sigmax); Fstep=(Fmax-Fmin)/(I-1);
x(1)=xmin; currentF=Fmin;
for i=2:I-1
    currentF=Fmin+(i-1)*Fstep;
    x(i)=Finv(currentF,myx,sigmax);
end
x(I)=xmax;

z=zeros(1,H);
Fmin=F(zmin,myz,sigmaz); Fmax=F(zmax,myz,sigmaz); Fstep=(Fmax-Fmin)/(H-1);
z(1)=zmin; currentF=Fmin;
for h=2:H-1
    currentF=Fmin+(h-1)*Fstep;
    z(h)=Finv(currentF,myz,sigmaz);
end

% for h=2:H-1
%     currentF=Fmin+(h-1)*Fstep;
%     z(h)=Finv(currentF,myz,sigmaz);
% end
z(H)=zmax;

dx=diff(x); dz=diff(z);